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ABSTRACT 



This thesis develops the concept of minimum-time (Bang-Bang) controllers 
and their application to missile control. Based on Pontryagin’s minimum 
principle, a minimum-time second order controller is derived. This controller is 
then applied to control of a vertically-launched surface-to-air missile. In the 
boost phase of missile flight, the minimum-time controller drives the missile 
body axis from vertical to a commanded angle in minimum time. In the terminal 
phase of the missile-target engagement, the minimum-time controller drives the 
time rate of change of the line of sight angle to zero in minimum time. The 
results obtained with the minimum time controller are compared with those 
obtained with the Proportional Navigation control algorithm, which is 
commonly used in tactical surface-to-air and air-to-air missile. 
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DISCLAIMER 

The reader is cautioned that computer programs developed in this research 
may not have been exercised for all cases of interest. While every effort has 
been made, within the time available, to ensure that the programs are free of 
computational and logical errors, they cannot be considered validated. Any 
application of these programs without additional verification is at the risk of the 
user. 
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I. INTRODUCTION 



In modem missile warfare, new technologies are producing faster, more 
accurate anti-ship cruise missiles (ASCM). They are able to fly at very low 
altitudes and perform a high speed 'pop-up' maneuver just before impact, in 
order to defeat the ship's missile defenses, as shown in Figure 1. In order to 
counter the advanced ASCM's, new control schemes must be developed to 




Figure 1. Anti-ship cruise missile terminal maneuver 



enable the ship's surface-to-air missiles to counter this threat. 

The aim of this thesis is to illustrate the concept of Bang-Bang control and 
how it applies to control of vertically-launched surface-to-air missiles. Bang- 
Bang control can be used in the missile boost phase in order to turn the missile 
downrange in minimum time, and also in the terminal phase, where the missile 
must react very quickly in order to intercept a maneuvering target. 

This work applies the minimum-time (Bang-Bang) control concept to 
vertically-launched surface-to-air missiles. In Chapter II, the concept of 
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control is defined and solved for a second order system. In Chapter HI, the 
equations of motion for a vertically-launched missile in the boost and terminal 
phases are developed. In Chapter IV, the boost simulations comparing missile 
trajectory versus maximum thrust vector control (TVC) angle are presented. 
In Chapter V, the terminal phase missile simulations are presented for a non- 
maneuvering target and also for a maneuvering target. 
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II. CONCEPT OF BANG-BANG CONTROL 



In this section the form of the optimal control for a particular class of 
systems will be determined by using Pontryagin’s Minimum Principle. 

A. PONTRYAGIN’S MINIMUM PRINCIPLE 

It shall be assumed that the state equations of the system are of the form 

x = Ax + Bu, (2-1) 

where A is an n by n array, B is an n by m array, and A and B may be 
explicitly dependent on the states and time. For this research it will be assumed 
that there is a single input; therefore, m is 1. It is specified that the admissible 
controls must satisfy the inequality constraints 

N. < u(t) < N+ (2-2) 



where N. and N+ are the known lower and upper bounds for the control input. 
It is desired to drive x from an initial state x(to) to to a desired final state x(tf), 
where to is the problem start time and tf is the problem end. The optimal control 
is that control which drives the state from its initial state to desired final state 
using the least "cost". The cost can be any desired measure of system 
performance. In the case of a minimum-time controller, the cost function can 
be represented as 



tf 

J = Jdt = tf-to 

^0 



(2-3) 
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Pontryagin's Minimum Principle states that the optimal control, u*, which 
minimizes the cost function, must minimize the Hamiltonian, which is defined as 



H(x(t), u(t), p(t), t)= l+pT(t)(Ax(t)+Bu(t)) (2-4) 



where p(t) is a Lagrange multiplier vector, and is arbitrary. p(t) can be defined 
as the 'costates' of the system, and p(t) can be written as 



P(t) = 



P2(0.' 



(2-5) 



The necessary condition for the optimal control u* to minimize the cost 
function J is 



H(x*(t),u*(t),p*(t),t) < H(x*(t),u(t),p*(t),t) (2-6) 



for all times between to and tfand for all admissible controls. The asterisks in 
equation (2-6) represent the optimum values, and thus. 



l+p*'T(t)Ax*(t)-i-p*'^(t)Bu*(t) >l-i-p*T(t)Ax(*t)+p*'^(t)Bu(t) (2-7) 



Therefore 



p*T(t)Bu*(t) < p*T(t)Bu(t) (2-8) 



for all admissible u(t) and for all times from ti„i,iai to tfi„ai. If u(t) is constrained to 
±N, the optimal control U*(t) is 



u*(t) = 




for P 2 > 0 
for P 2 < 0 



= -Nsign(p2) 



(2-9) 
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This indicates that the time-optimal control is 'Bang-Bang'; that is, the 
optimum control switches between its maximum positive and negative values. 
There are also three ideas which deserve mention: 

(1) (Existence) 

If all the eigenvalues of A have nonpositive real parts, then an optimal 
control exists that transfers any initial state xq to the origin. If there 
are positive eigenvalues (i.e. unstable roots), there may be some 
region of the state space where the system is uncontrollable. 

(2) (Uniqueness) 

If an optimal control exists, then it is unique. 

(3) (Number of Switchings) 

If the eigenvalues of A are all real, and a (unique) time-optimal 
control exists, then the control can switch at most (n-1) times. 

Thus, an nth-order system having all real, non-positive eigenvalues has a unique 
time-optimal control that switches at most n-1 times (not counting switching off 
at time tfmai)- For complex conjugate eigenvalues, more than n-1 switchings 
may be required. [1] 

B. APPLICATION TO A SECOND ORDER SYSTEM 

In order to illustrate the minimum principle, consider the second order 
differential equation 

x = u, (2-10) 



or in state space. 
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( 2 - 11 ) 
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From the minimum principle, the time-optimal control for this system is +N. 
Thus the segment of optimal trajectories can be found by integrating the state 
equations (2-11), with u=±N, from time to to t. This yields the results 

X2(t) = +Nt + Ci (2-12) 



Xi(t) = ±N — -HCit + C2 
2 



(2-13) 



where 



Cj = +Nto + X2(to) 



(2-14) 



and 



C2 --i-N-^-Cito -i-Xi(to). 



(2-15) 



Cl and C 2 are functions of the initial conditions and can be treated as constants. 
Time can be eliminated from equations (2-12) and (2-13) by squaring the first 
equation, dividing the result by 2N, and comparing the result with equation 
(2-13) to obtain 

Xi(t) = — X 2 ^(t)+C 3 forU = +N (2-16) 

2^ 



Xl(0 = -^X2^(t)-HC4 



forU = -N 



(2-17) 



where 
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(2-18) 



Cs 




C4 



= C. 



(^1 

2N 



(2-19) 



C3 and C4 are functions of the initial conditions and may also be treated as 
constants. Equations (2-18) and (2-19) each define a family of parabolas, as 
shown in Figures 2 and 3. 

Analysis of Figures 2 and 3 reveal the controls corresponding to the 
following situations: 

1. u*=+l implies that the initial state Xq lies on segment AO at time to. 

2. u*=-l implies that the initial state xq lies on segment BO at time to. 

This defmes the optimal control to reach the origin. In order to reach these 
zero trajectories, the opposite value of control is required. It can be seen from 
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Figure 4 that the two segments AO and BO form a single curve. This curve is 
known as the switching curve. All states that fall above the switching curves 
will result in the control input U=-N, and those states which fall below the 
switching curves will result in the control input U=+N. This curve can be 
represented as 

Xi(t) = -^X2(t)|x2(t)| (2-20) 



Thus, the control law can be obtained by moving all the terms of equation 
(2-20) to one side: 



U(t) = -Nsignfxi(t)-H^X2(t)|x2(t)| 



( 2 - 21 ) 
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Examples of the use of this switching law is shown in Figure 5. From any 
arbitrary initial state, the trajectory will follow the parabola that passes through 
the start point until the trajectory reaches the optimal switching curve. Then the 
control switches signs, and follows the optimum parabolic path to the origin. 

Equation (2-21) is the mathematical representation of the Bang-Bang 
control law. This control law can be used in numerous applications where it is 
desired to transfer a state from initial conditions to desired final conditions in 
minimum time. Two examples that are of interest in this research are: 

1) drive a state, such as a position, to zero in minimum time. The control 
law for this portion of the missile flight is of the form 

U(t) = -Nsign|^Xi(t) + ^X2(t)|x2(t)|j, (2-22) 
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where Xi(t) is the position at time t, X 2 (t) is the velocity at time t, and N is the 
maximum admissible force that can be used to accelerate or decelerate the 
state; and 

2) drive the time rate of change of a state to zero. In this case it is 
desired to keep a state constant, as opposed to driving it to zero. The form of 
the controller is also given by equation (2-22), but X[ in this case will represent 
velocity and X 2 will represent acceleration, with N being the maximum 
admissible control with which to change the states. 

In the next chapter, the application of the second order Bang-Bang 
controller to boost phase and terminal phase missile control will be examined. 
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III. DEVELOPMENT OF EQUATIONS OF MOTION FOR A 
VERTICALLY-LAUNCHED MISSILE 



The equations of motion for a missile can be derived from Newton's second 
law of motion, which states that the summation of all external forces acting on a 
body must be equal to the time rate of change of its momentum, and the 
summation of the external moments acting on a body must be equal to the time 
rate of change of its moment of momentum (angular momentum). The time 
rates of change can be expressed by two vector equations: 



where i indicates i‘*» cartesian coordinate of the vector with respect to inertial 
space. By definition, H is the angular momentum, or moment of momentum, of 
a revolving body. [2] 

There are several assumptions that are made in order to simplify the 
problem: 

1. The mass of the missile remains constant. 

2. The missile airframe is a rigid body. 

3. The earth is an inertial reference, and unless otherwise stated, the 
atmosphere is fixed with respect to earth. 






(3-1) 



and 




(3-2) 
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In order to completely describe the motion of a missile in three dimensions, 
a total of six nonlinear differential equations must be solved. However, certain 
assumptions will be made in order to reduce the number of equations needed to 
adequately describe the motion of the missile: 

1. The X axis will be assumed to be downrange from the launch 
platform to the target, the Y axis will be crossrange, and the Z axis will be 
altitude. 

2. The type of missile studied here is known as a skid-to-turn missile, 
as this missile uses direct side force to turn. Thus, there will be no coupling 
terms between the pitch and yaw axes, and one set of equations will describe 
the motion in each of the Y and Z planes. 

3. It will be assumed that the missile is roll stabilized, i.e., that there 
will be no turning moments about the X axis. 

4. Since there are no cross-coupling terms between the pitch and yaw 
axes, the simulation can be broken into two problems, one where the pitch 
angle is held constant, and the other where the yaw angle is held constant. In 
this discussion the yaw angle will be assumed to be constant at zero. Using 
these assumptions, a 3 degree-of-freedom model can be developed. The force 
diagram of a 3 DOF missile in flight is shown in Figure 6. The governing 
equations of motion of this missile are: 

2F^=mx= -FiiftSiny^ -i- Tcos(0^ + 5) - F^ragCosYm (3-3) 
2 Fz = mz = - mg + Fiift cosy^, -i- T sin(0jn -h 5) - sin (3-4) 

^ M(,g — Icg0 ~ ~ FjjfjLpg cos'Yjn ~ Fjjj. 3 gLpg sin Yjj, — TL(g sin 5 t3-5') 
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where: 

m missile mass (Kg) 

g gravity (9.81 m/s^) 

Fjift aerodynamic lift force (N) 

Fdrag aerodynamic drag force (N) 

Tm angle of velocity vector from reference 
0m angle of missile head from reference 

5 TVC deflection angle 

leg moment of inertia of the missile about the center of gravity 
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Lpg distance from aerodynamic center of pressure to center of 
gravity 

Ltg distance from TVC control to center of gravity 
The aerodynamic force acting on a lifting body, in this case the wings and 
control surfaces of the missile, may be resolved into two components, which 
are 

Fiif.=CLfv^S (3-6) 

and 



P —C' Q 

^drag ^ 



(3-7) 



where (p/2)Vm^ is the dynamic pressure in N/m^. 

The lift and drag coefficients are functions of angle of attack (a), and Mach 
number. Actual values for Cl and Cd must be measured, as they are highly 
dependent on the airfoil geometry. However, Cl and Cd may be modeled 
satisfactorily using the following approximations: 



Cl =0.1a 



(3-8) 



Cd = 



2T 



pVmaxS 



(l-i-0.2a^) 



(3-9) 



These coefficients are in reasonable agreement with reference [2] at low 
angles of attack (a<10°). The only time the missile can expect to exceed 10° 
angle of attack is during the tip-over phase, when velocity, and therefore lift and 
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drag forces, are small compared to the thrust force. This is considered an 
appropriate approximation for the scope of the model under discussion. 



A. BOOST PHASE CONTROLLER 

The primary objective of the missile control system in the boost phase of 
flight is to tip the missile over and point the missile in the general direction of the 
target in minimum time. Therefore, 0, the angle that the missile makes with the 
horizontal reference, is a logical choice for a control variable. If the force that 
creates the moment that tips the missile over can be controlled, a general 
control algorithm can be formulated. This is shown graphically in Figure 6. 

At launch, the missile is assumed to be in a vertical position, with the 
velocity vector pointing straight up, i.e. 0 = 90° and y = 90°. In the initial boost 
phase, it is assumed that the first two terms of equation (3-5) are negligible as 
compared to the third. Therefore, equation (3-5) may be simplified to 



0 = 



-TL(.g sin 5 



‘eg 



(3-10) 



In state space, this can be represented as 



’ 0 ' 




'0 
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■0‘ 


- 1 - 


'O' 


^-TLpgSin5^ 
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o_ 
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V eg y 



(3-11) 



The only term in equation (3-11) that can be used to control 0 is 5, the TVC 
deflection angle. A simple and effective method for controlling missile tip over 
is to employ the Bang-Bang controller, developed in Chapter II, to drive 0 to 0 
(horizontal flight) in minimum time, as shown in Figure 7. Thus, the control law 
is 
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f 1 

u(t) = -Nsign|^0(t) + — 




( 3 - 12 ) 



where N is the maximum torque available to rotate the missile and is 
determined from the relation 



N = -IbLsin(5„3,) 



I 



eg 



( 3 - 13 ) 



and 6inax is the maximum thrust vector control (TVC) angle available. 

This indicates that for minimum tip-over time, 8max=90°, However, this 
may not be desirable since all of the thrust will go to rotate the missile, while 
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none of the thrust will accelerate the missile downrange. This amounts to the 
missile 'pin- wheeling' in the air. Therefore, TVC angles of less then 90° should 
be used. Several TVC angles will be examined to show the tradeoffs involved. 

When the missile achieves horizontal orientation (i.e. 0=0°), the control law 
will continue to switch from -i-N to -N in an attempt to keep the missile at 
exactly 0°. This phenomenon is known as 'chattering'. At this point it is no 
longer desirable to use maximum control effort to correct for small changes in 
0. A satisfactory method for removing the chatter from a Bang-Bang 
controller is to employ linear control in a small region, of width e, about the 
desired value, as shown in Figure 8. 

This adapted Bang-Bang controller can be expressed mathematically as 



u(t) = 



■Nsignre(t)+^e(o|e(t)| 



e(t)+^6(t)|e(t) 



>e 



-^fe(t)+^e(o|e(t)|l -e<e(t)+^e(t)|e(t)|<e. o-m) 



J 

4N- 

-e 


+U(t) 

/ 1 ^ i 4 /^\ 


k 

1 


1 

-N 



Figure 8. Bang-Bang Control Law with Linear Zone 
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B. TERMINAL PHASE 

1. Missile-Target Geometry 

In this section, it will be assumed that the vertically launched surface- 
to-air missile has completed tip-over, and has settled out to a constant speed 
and altitude. The start point for the terminal engagement of the target begins 
with the missile seeker acquiring the target. The geometry of the missile-target 
engagement is shown in Figure 9. 



In describing the geometry of the missile-target encounter, several 
parameters must be calculated. The angle the missile's velocity vector makes 
with the inertial reference frame is defined as y: 




yt-y. 



m 



Figure 9. Missile-target geometry 




(3-15a) 
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The magnitude of the missile velocity vector Vm is given by the relation: 



V = +v^ 



(3-15b) 



where Vym is the vertical component of the missile velocity vector and Vxm is 
the horizontal component of the missile velocity vector, measured in inertial 
coordinates. 

Similarly, the magnitude of the target velocity vector and and angle it 
makes with the inertial reference is given by 



In computing the variables used to calculate the line of sight and its 
derivatives, it is convenient to describe the relative position, velocity, and 
acceleration between the target and the missile: 




(3- 16a) 



and 




(3-16b) 



X = XfXm 



(3-17a) 



Y = Yt-Ym 



(3-17b) 



Vx = Vxt-Vxm 



(3-17c) 



Vy — Vyf Vym 



(3-17d) 
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Ax = Axt - Axm 



(3-17e) 



Ay — Ayt - Aym. (3"17f} 

The line of sight between missile and target is given by O, and is 
defined as: 

(3-18) 

where X and Y are the relative cartesian positions as defined by equations (3- 
17a) and (3- 17b). The first and second time derivatives of a, denoted <j and a, 
are computed using the relations 

^ (xVy-YV,) 

(3-19) 

and 

2(XV,-hYV,)(XV,-YV,) 

R'* 

In an actual missile-target engagement, it is very unlikely that the 
target line of sight rate and acceleration will be measured accurately enough in 
order to use d and d directly in the guidance algorithm. An alternative 

. /V 

method is to calculate estimates of d and d, denoted by d and d , based on 
measured values of a. This is accomplished through the use of a Luenberger 
Observer, as shown by the signal flow graph in Figure 10. 



a = tan' 



-ifY' 

lx. 
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Figure 10. Luenberger Observer for estimating a and a 
This observer can be described in state-space form by the relation 
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( 3 - 21 ) 



where R is the range from the missile to the target, R is the range rate, and ki, 
k 2 , and k 3 are gains to be determined by the designer in order for the estimated 
states to closely follow the true states. 

2. Missile Guidance Laws 

The objective of tactical missile guidance is to keep the line of sight 
angle, a, between missile and target constant. This is desirable because if the 
line of sight remains constant, while the range from missile to target decreases, 
an intercept will occur, as shown in Figure 11. Two missile guidance laws will 
be presented: classical proportional navigation and minimum time (Bang-Bang) 
control. 
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a. Proportional Navigation Control 

The method most often used in controlling homing missiles is 
known as proportional navigation, where the magnitude of the transverse 
control force is proportional to the rate of change of the line of sight (6). 
Lukenbill [3] conducted research into proportional navigation control and his 
derivation of the proportional navigation controller is given below. 

Figure 12 depicts the basic proportional navigation scheme. 
Assuming that the seeker head of the missile follows the target, the transverse 
acceleration perpendicular to the line of sight will equal the acceleration of the 
R vector in that direction. Mathematically, the acceleration of R is 

Ar =(R+coxcoxR)iR +(2coxR + cbxR)ie (3-22) 



where 

R 



missile/target line of sight vector 
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R 



closing rate along R 
R = acceleration along R 

CO = angular rate of change of R in inertial space 

Ar = overall acceleration of R. 

At this point, a missile acceleration, A„, equal to the target 
acceleration. A,, will make the line of sight parallel to its original direction. As 
long as R remains along R (co=0) a missile/target intercept is assured. So, the 
transverse acceleration command is 

At - Ajj, =coxR + 2((oxR) (3-23) 
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Figure 12. Vectorial Proportional Navigation Scheme 



Assuming the line of sight rate is equal to the angular rate of change of R in 
inertial space, equation (3-23) now becomes 



Af-An, = Ra+2Ra 



(3-24) 



where d equal to co, and d is equal to cb. 

In the classical proportional navigation scheme, the missile course 
is one in which the rate of change of the missile heading is directly proportional 
to the rate of rotation of the line of sight vector from the missile to the target. 
As a result, this course change is intended to counteract the rotation of the line 
of sight, thus returning to a constant bearing course. The movement of the 
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missile and target cause the line of sight to rotate, resulting in a differential 
displacement between the missile and the target perpendicular to the range line. 
Figure 13 depicts this geometry. The proportional navigation guidance law 
attempts to generate an acceleration command, Ac, perpendicular to the line of 
sight 

Assume a gyro stabilized seeker head, as in the Sidewinder 
missile. If there is no torque applied to the gyro, the seeker will not rotate. 
Assuming the seeker tracks the target, the gyro angle will follow the line of 
sight. Applying the equation of motion for a gyro stabilized seeker 

L=I(oQ (3-25) 



where 


L 


applied torque 




CO = 


spin angular velocity 




I 


moment of inertia of the gyroscope 




Q 


rate of precession of the gyroscope. 



Applying this to the case when the seeker head tracks the target, Q is then 
replaced by the rate at which the gyro is torqued in space. This is simply a, 
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Vm(t) \ 

/ ^X^Vm(t+dt) 



Figure 14. Missile Acceleration Relationship 
which is the line of sight rate. Thus, equation (3-25) becomes 

L = I(od. (3-26) 

This torque is in turn applied to the control surface of the missile leading to the 
relationship 

= kL = kicoa (3-27) 

where k is a constant of proportionality. Referring to Figure 14, a relationship is 
determined for Am in terms of the rate of change of the missile flight path angle, 
Ym- Given the missile velocity vector at some point in time, V„,(t), and suppose 
the missile undergoes an acceleration, A„, during an interval of time, dt. The 
velocity vector is then displaced and is represented by the vector V„(t+dt). The 
angle the vector is traversed is simply dYm» the differential missile flight path 
angle. For small angles (which are guaranteed by making dt small) the 
following relationship is obvious: 



Amdt=V^dYm- 



(3-28) 
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Dividing equation (3-28) by the time interval, dt, the missile acceleration is 
defined as 



A =v = V Y 

’mim. 



(3-29) 



Combining equations (3-27) and (3-29) 



Vn,Tm = kIo^. 



(3-30) 



Dividing through by Vm, the proportional navigation law becomes 



f 



Tm = 



kIcoY 






J 



(3-31) 



or 

tm = (3-32) 

Equation (3-32) represents the classical proportional navigation equation where 
Ynj = rate of change of the missile heading 

a = rate of change of the line of sight 

N = proportional navigation ratio. 

The navigation ratio determines the sensitivity of the missile 
system. A high navigation ratio will lead to rather high gains resulting in large 
missile commands for small changes in the line of sight rate. On the other hand, 
small values for N will lead to small missile commands for a given d. Larger 
navigation ratios are preferred for head on engagements and smaller ones are 
preferred for tail chase cases. For this research the navigation ratio is taken to 
be four. 
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In tactical radar homing missiles using proportional navigation 
guidance, the seeker provides an effective measurement of the line-of-sight 
rate, and a Doppler radar provides closing velocity information. [4] 

b. Bang-Bang Control Minimizing Line-of-Sight Rate 

In this control law, the goal is the same as in proportional 
navigation control: to drive the line-of-sight rate to zero, or in other words, to 
keep the line of sight constant. 

Mathematically, this control law can be stated as. 



u(t) = -N sign 



a|a 



(3-33) 



In order to alleviate chattering, a small linear region will be introduced, similar 
to the boost phase controller: 



u(t) = 



-Nsign t) -h ^ a( t )|a( t)| j 



ct(t)-f-^a(t)|a(t)| 



>e 



-Y[a(t)+^d(t)|a(t)|j -e<a(t)+^a(t)|d(t)|<e.(3-34) 



In the next chapter, several simulations will be presented in order 
to illustrate the effectiveness of the Bang-Bang control algorithms as compared 
to Proportional Navigation, in both the boost phase and in the terminal phase of 
the missile flight. 

3. Missile Dynamics 

The signal flow graphs for the missile dynamics are shown in Figure 
15. 
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In state-space, these signal flow graphs can be represented as 
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IV.VERTICAL LAUNCH MISSILE BOOST PHASE 

SIMULATION 



A. OVERVIEW 

In simulating the boost phase of the vertically launched missile, several 
assumptions can be made to reduce the complexity of the simulation to a level 
that won’t detract from the concepts being presented: 

1 . The missile is limited to Mach 4. 

2. The missile is limited to 30 g's of acceleration in the transverse and 
tangential directions combined. 

3. The speed of sound in air is constant at 340 m/s. 

4. The density of air is constant at 1 Kg/m^. 

5. The missile mass remains constant. 

6. The lift and drag forces experienced by the missile can be modeled by 
equations (3-6) through (3-9). 

B. MISSILE PARAMETERS 

At launch, the following parameters define the state of the missile: 

The initial missile parameters are: 



Lm=4.2 m 


missile length 


M=225 Kg 


. missile mass 


Ltg=2.1 m 


length from center of gravity to tail 


Lpg=0.15 m 


length from center of gravity to center 




of pressure 


S=0.99 m2 


missile reference lift surface 
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T=66800 N 

I,g=2000 kg m2 
p=l. kg/m3 
V„,^=1460m/s 



missile thrust 

missile moment of inertia 

density of air 

max missile velocity 



Additionally, the launch point will correspond to inertial reference: 
x(0)=0 m y(0)=0 m 

Vx(0)=0 m/s Vy(0)=0 m/s 

C. CANNISTER EGRESS 

It is not possible, nor is it desirable, to have the missile commence its tip- 
over maneuver until it has egressed from the launch cannister and achieved 
sufficient distance from the launch platform to minimize hazards to personnel. 
Thus, the TVC actuators will have a zero degree deflection for the first 0.7 
seconds after launch. This will ensure sufficient missile altitude before tip-over 
begins. At time=0.7 seconds the missile position and velocity are: 
x(0.7)=0m y(0.7)=71.8m 

V^(0.7)=0 m/s Vy(0.7)=204.5 m/s 

D. BOOST SIMULATION AND RESULTS 

The simulation was run with four values of the maximum thrust vector 
control (TVC) angle: 30°, 45°, 60°, and 90°. This was done to illustrate what 
effect a greater TVC angle will have on the missile trajectory and velocity. 
Figure 16 shows the missile trajectories for the four cases. 
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Figure 16. Missile Boost Phase Trajectories 



Analysis of the trajectories shows that as the maximum TVC angle is 
increased, the missile will achieve horizontal flight at lower altitudes than 
smaller TVC angles. However, at angles greater than 45°, there is a significant 
reduction in the downrange distance. This may be crucial if the missile is to be 
fired against a target that is close to the ship at launch. At TVC angles less 
than 30°, the thrust is directed toward achieving missile altitude to a larger 
extent than having the missile achieve horizontal flight in minimum time. This 
may be advantageous if the missile if fired against a high-altitude target, but not 
against a sea- skimming missile. Another consequence of having a large TVC 
angle is that much of the thrust is spent rotating the missile to horizontal at the 
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expense of missile kinetic energy, which results in a longer time for the missile 
to reach maximum velocity. This is shown graphically in Figure 17. 




Figure 17. Missile Boost Phase Velocity Magnitude Profiles 

At the two highest TVC angles (60° and 90°), the first several seconds are 
spent turning the missile over while not attaining much forward velocity. As a 
result of a higher turning rate, higher angles of attack are generated, as shown 
in Figure 18. This is not desirable because at higher angles of attack, the drag 
forces that act to slow the missile down increase. Thus, it is advantageous 
from an aerodynamic point of view to keep the angle of attack as small as 
possible in order to minimize the drag force acting on the missile. 

Figure 19 shows the missile pitch angle as a function of time. It is of 
interest to note that there is not a great difference between total turning times 
between the 30° TVC angle and the 90° TVC angle. This is due to the counter- 
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acting lift and drag forces encountered by the missile, during high TVC angles, 
that are caused by high angles of attack. This indicates that increasing the 
maximum TVC angle beyond 45° may achieve a somewhat faster turning time, 
but at the expense of a great loss of kinetic energy and therefore velocity. This 
is of peak importance, for without sufficient velocity, a missile has little or no 
chance of successfully intercepting an incoming target, particularly a 




Figure 18. Missile Boost Phase Angle of Attack Profiles 
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Figure 19. Missile Boost Phase Pitch Angle ProHles 
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V. TERMINAL PHASE SIMULATIONS 



In order to compare the effectiveness of minimum-time 'Bang-Bang' control 
as compared to proportional navigation control, several missile-target 
geometries were simulated. Each geometry was run with four control system 
configurations: 

1. proportional navigation using analytical values of line-of-sight rate and 
acceleration, which shall be referred to as 'analytical proportional navigation’; 

2. proportional navigation using estimates of line-of-sight rate and 
acceleration, which shall be referred to as 'estimated proportional navigation’; 

3. Bang-Bang control using analytical values of line-of-sight rate and 
acceleration, which shall be referred to as 'analytical Bang-Bang’; 

4. Bang-Bang control using estimates of line-of-sight rate and acceleration, 
which shall be referred to as 'estimated Bang-Bang'; 

A. CASE ONE 

In the first set of simulations, the target will be headed directly toward the 
ship, flying a straight course with no maneuvers. The initial conditions are: 



Missile: 



x=0; 

Vx=680 m/s (Mach 2) 



z= 1000 m; 
Vz=0; 



Target: 



x=7000 m; 



z=30m 



Vx=-850 m/s (Mach 2.5); Vz=0; 
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Figure 20 shows the missile-target geometries using true line-of-sight rate 
and acceleration. Both proportional navigation and Bang-Bang control was able 
to hit the target. Note that the trajectory using Bang-Bang control flew a 
straighter course to the intercept point than the trajectory using proportional 
navigation control. This indicates that Bang-Bang control is able to guide the 
missile onto an intercept trajectory faster than the proportional-navigation 
controller, which is continuously changing the missile flight path in order to 
intercept the target. 

Figures 21 and 22 show control force as a function of time for the 
proportional navigation controller and the Bang-Bang controller, respectively. 
The proportional navigation controller required continuous control in order to 




Figure 20. Case One: Missile/Target trajectories using true line 

of sight rate 
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correct the missile trajectory enough to enable the missile to intercept the ; 
target. It does not leave much reserve control to counter a maneuver away , 

from the missile by the target. Conversely, the Bang-Bang controller applied ^ 

maximum control effort until the line of sight rate became zero, then shut off. 

This allows the missile greater flexibility in case the target maneuvers in any 
direction. If the linear zone were not incorporated into the Bang-Bang 
controller, however, the controller would have chattered between maximum 



positive and negative control force instead of shutting off. This would have 
decreased the missile velocity unnecessarily just when the missile needs 
maximum kinetic energy for maneuverability. 




Figure 21. Case One: Control force versus time using proportional 
navigation control with true line of sight rate 



Figure 23 shows the trajectories using estimates of line of sight rate and 
acceleration. It is apparent here that the proportional navigation controller was 
not able to react fast enough in order to hit the target. This is due to the delay 
encountered while the observer filter was not matched up with the analytical 
values, as shown in Figure 24. By the time the estimated value of line of sight 
rate and acceleration matched the analytical values, the target was at a point 
where proportional navigation could not generate enough of a control signal fast 
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enough to enable target intercept. This illustrates the importance of being able 
to detect, virtually instantaneously, changes in target velocity and direction.The 
Bang-Bang controller also experienced the same delay before the proper 
control was applied, as shown in Figure 25, but since it immediately applied 
maximum control, the missile was still able to intercept the target. 
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Figure 23. Case One: Missile/Target trajectories using estimated 

values for line of sight rate 
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Figure 24. Case One; Control force versus time using proportional 
navigation control with estimated values for line of sight rate 




Figure 25. Case One: Control force versus time using Bang-Bang 
control with estimated values for line of sight rate 

Figure 26 compares the trajectories of Bang-Bang control using analytical 
and estimated values of line of sight rate and acceleration. Because the 
estimated values took approximately 0.25 seconds to match the analytical 
values, there is a period at the start of the simulation where the estimated 
trajectory has positive control applied to it, while the analytic trajectory has 
negative control applied. After the analytic and estimated values matched up, 
they both were able to drive the line of sight rate to zero in plenty of time to 
intercept the target. 
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Figure 26. Case One: Missile/Target trajectories using Bang-Bang 
control, comparing true and estimated values for line of sight rate 



B. CASE TWO 

In the second set of simulations the target will be headed directly toward 
the ship, but will initiate a 12 g vertical maneuver at problem start. The initial 



conditions are: 




Missile: 




x=0; 


z=10CX) m 


Vx=680 m/s (Mach 2) 


Vz=0; 


Target: 




x=7000 m 


z=30m 


Vx=-850 m/s (Mach 2.5); 


Vz=0. 
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Figure 27 shows the simulation results of the analytical proportional 
navigation and analytical Bang-Bang controllers. The Bang-Bang controlled 
missile was able to adjust to the 12 g vertical maneuver of the target in order to 
achieve intercept. The proportional navigation controlled missile was not able 
to keep up with the high g maneuver of the target, and missed the target by 
about 80 meters. 

Figure 28 shows the control force as a function of time for the analytical 
proportional navigation controlled missile. Note that the controller did not 
saturate until well into the simulation. This indicates that if more control had 
been applied sooner, the missile could have intercepted the target. This is 
illustrated by Figure 29, which shows the control of the analytical Bang-Bang 




Figure 27. Case Two: Missile/Target trajectories using true line of 

sight rate 
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Figure 28. Case Two: Control force versus time using Proportional 
Navigation control with true line of sight rate 




Figure 29. Case Two: Control force versus time using Bang-Bang 
control with true line of sight rate 



controlled missile. Neither controller went to zero as the Bang-Bang controllers 
of Case One did. This is because the target is maneuvering, constantly 
changing the line of sight angle. 

The simulations using estimated values for line of sight rate and 
acceleration are shown in Figure 30. As was the case with the analytical 
expressions, the Bang-Bang controller using estimates was still able to achieve 
intercept, while the proportional navigation controller missed the target. Note 
that the miss distance for the proportional navigation controller in Figure 30 is 
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less than the miss distance achieved with the analytical expressions. This is 
because the target turned toward the missile, and the delay caused by the 
observer prevented the missile from beginning to dive after the target as 
quickly, thus leaving less of an altitude difference when control is applied. If the 
target were to turn away from the missile, the estimates of line of sight rate and 
acceleration would be worse than those obtained using analytical values. 

Figures 31 and 32 are the control versus time graphs for the estimated 
proportional navigation controller and the estimated Bang-Bang controller. The 
delay encountered by both of the controllers is about the same as those 
encountered in Case One. 




Horizontal Range (meters) 



Figure 30. Case Two: Missile/Target trajectories using estimated 

values for line of sight rate 
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Figure 31. Case Two: Control force versus time using proportional 
navigation control with estimated values for line of sight rate 




Figure 32. Case Two: Control force versus time using Bang-Bang 
control with estimated values for line of sight rate 
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VI. CONCLUSIONS AND RECOMMENDATIONS 



This work has shown that the minimum-time (Bang-Bang) controller is an 
effective algorithm for missile control, in both the boost and terminal phase of 
tactical missile flight. It is particularly effective when the target has a speed 
advantage over the missile or when the target is maneuvering. 

As is the case with virtually all tactical missiles, the accuracy of the target's 
measured position, velocity, and acceleration vectors are paramount in 
accurately predicting the parameters with which the missile is controlled with a 
Bang-Bang minimum time controller. 

Areas for future study will include: 

1. Simulating the entire missile flight from launch through intercept. The 
focus of these studies will be to analyze the effects of missile velocity, tip-over 
altitude, and target speed advantage on missile performance. 

2. Developing more complex models in order to better understand what 
the effects of a Bang-Bang controller will have on a particular missile system. 
Such studies will use empirical aerodynamic and physical data for a particular 
missile, such as SM-II Block FV or Vertically-Launched Sea Sparrow, and will 
entail developing a 5 or 6 DOF computer model. 

3. Investigating the effects of noise on measurement of line of sight and 
its time derivatives, to the extent of the effect noise has on a Bang-Bang 
controller. 

4. Investigating the effects of plant modelling errors on system sensitivity 
when using Bang-Bang controllers. 
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APPENDIX 1- BOOST PHASE MATLAB PROGRAM 



% VLS BOOST PHASE MISSILE SIMULATION 
% T B MULL 

% REVISED 22 APRIL 1992 
clear;clg; 

% Initial Missile Parameters 
Delta=[30 45 60 90]; 
for j= 1 : length(Delta) 



L_m=3.6; 


% 


L_m= missile length in meters 


M=227; 


% 


M= missile mass in kilograms 


L_gt=1.8; 


% 


L_gt= length from center of gravity to tail 


L_gp=0.15; 


% 


L_gp= length from center of gravity to center 




% 


of pressure 


S=0.99; 


% 


S is lift surface in square meters 


T=66800; 


% 


T= missile thrust in newtons 


Icg=2000; 


% 


leg is the missile moment of inertia (kg m'^2) 


rho=1.00; 


% 


rho is the density of air in kg/m'^3 


Vmax=340*4; 


% 


Vmax is the max missile velocity in m/s 


dt=.02; 


% dt is simulation time step size 


delta_max=Delta(j); 







% 

% Initial Conditions 



tfinal=10.00; % Final simulation time. 
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kmax=tfinal/dt; 


% kmax is number of time increments 


time=0.0; 


% initial time 


x=0;xdot=0;xddot=0; 


% initial missile states (in x) 


y=0;ydot=0;yddot=0; 


% initial missile states (in y) 

% (these initial y states correspond to missile 
% position and velocity after clearing the 
% cannister) 


gamma=90.0; 


% gamma is inertial to velocity vector 


gammadot=0; 


% in degrees (will convert to radians in loop) 


theta=90.0; 


% theta is inertial to missile head 


thetadot=0; 




X_theta= [thetadot; the ta] ; 




X_x=[xdot;x]; 




X_y=[ydot;y]; 




A=[0 0;1 0] ; 


% A and B matrices for determining velocity and 
% position of x,y, and theta base on acceleration 
% commands 


B=[1;0]; 




[phi,del]=c2d(A,B,dt); 


% converts A,B to discrete for simulation 


alpha=theta-gamma; 


% alpha is Angle-of-attack 


dtr=pi/180; 


% converts degrees to radians 



V_m(l)=sqrt(X_x(l,l)'^2+X_y(l,l)'^2); % V_m is missile velocity magnitude 
F_drag=0; 

F_lift=0; 

Cd=0; 
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Q Q p 5tC 5fJ JfJ * * H! 5|5 H! * N« JfJ * * JfJ * HS H« JfJ JfS * JfS HS 

for i=l:kmax-l 



%****algorithm for determining rudder deflection angle****************** 

% test for cannister egress(time<0.6) 
if time(i)<0.7 
delta(i)=0; 
else 

Theta=X_theta(2,i); 

Thetadot=X_theta(l ,i); 

N=T*L_gt/Icg*sin(delta_max*dtr); 

test(i)=Theta+(Thetadot)*abs(Thetadot)/(2*N); 

eps=0.5; 



if abs(test(i))>eps 
delta(i)=delta_max*sign(test(i)); 
else 

delta(i)=delta_max/eps*test(i); 

end; 

end; 






F_Uft(:,i)=0.5*rho* V_m(i)^2*S*0. 1 *alpha(i); 

Cd(i)=(2*T/(M*rho*Vmax^2*S)); 

F_drag(:,i)=Cd(i)*0.5*rho*V_m(i)'^2*S; 
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%F_drag(:,i)=0; 

xddot(i)=-F_lift(i)/M*sin(gamma(i)*dtr)+T/M*cos(gamma(i)*dtr+delta(i)*dtr)... 

-F_drag(:,i)*cos(gamma(i)*dtr); 

yddot(i)=F_lift(i)/M*cos(gamma(i)*dtr)+T/M*sin(gamma(i)*dtr+... 

delta(i)*dtr)-F_drag(:,i)*sin(gainma(i)*dtr); 

thetaddot(i)=-F_lift(i)*L_gp/Icg*cos(alpha(i)*dtr)-T*L_gt*sin(delta(i)*dtr)/Icg; 

% 

X_theta(: ,i+l )=phi*X_theta(: ,i)4del*thetaddot(i); 
X_x(;,i+l)=phi*X_x(:,i)+del*xddot(i); 

X_y(:,i+l)=phi*X_y(:,i)+del*yddot(i); 

% - 

gamma(i+l)=atan2(X_y(l ,i+l),X_x(l ,i+l ))/dtr; 

% X_y=[ydot;y]; 

alpha(i+l)=X_theta(2,i+l)-gamma(i+l); 
delta(i+l )=delta(i); 

V_m(i+1 )=sqit(X_x(l ,i+l )'^2+X_y(l ,i+l T2)\ 

time(i+ 1 )=time(i)+dt; 

end; % END OF LOOP 

f*tf sic sic sic sic sic sic sIC sIC sIC sif sif sif sif sIc sIc sIc sIc *2^ 

•T* •T* *T* •T* •T* 5(5 5(5 5(5 ^5 5(5 

eval(['X_x',num2str(j),'=X_x;’]); 

eval([’X_y',num2str(j),’=X_y;’]); 

eval(['V_m',num2str(j),’=V_m;’]); 

eval(['alpha',num2str(j),'=alpha;']); 

eval(['X_theta’,num2str(j),'=X_theta;’]); 

delta_max 
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vm=V_m(i+l) 

X=X_x(2,i+l) 

Y=X_y(2,i+l) 

Alpha=alpha(i+1 ) 

Theta=X_theta(2,i+l ) 
end; 

Boostlplot2a 
save bocstdata 

clg; 

z=length(X_x(2,:)); 
for k=l :length(Delta) 
eval([’X_x=X_x',num2str(k), 
eval([’X_y=X_y',num2str(k), 
axis([0 9000 0 2000]); 
plot(X_x(2,:),X_y(2,:));grid; 
hold on 
ifk==l 

plot(X_x(2,l:2/dt:z),X_y(2,l:2/dt:z),’o’) 
plot(X_x(2,z),X_y(2,z),'o’) 
elseif k==2 

plot(X_x(2, 1 ;2/dt:z),X_y(2, 1 :2/dt: z),'o’) 
plot(X_x(2,z).X_y(2.z),’o’) 
else 

plot(X_x(2, 1 :2/dt:z),X_y(2, 1 :2/dt; z),’o') 
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plot(X_x(2,z),X_y(2,z),’o') 
end; 
end; 
hold off 

title(’Missile Boost Phase- Bang-Bang control'); 
xlabel('Horizontal Range (m)'),ylabel(’ Altitude (m)'); 

%text(0.6,0.8,’ 30® TVC angle', 'sc') 

pause;clg; 

axis('normal') 

% 

for k=l ;length(Delta) 
eval(['V_m=V_m’,num2str(k), ';’]); 
plot(time,V_m); 
hold on; 
end; 
hold off 

title('V_m versus time');grid; 
xlabelCTime (sec)’),ylabel('Velocity (m/s)’); 

% 

pause;clg; 

for k=length(Delta):-l:l 
eval(['alpha=alpha',num2str(k), ’;’]); 
plot(time, alpha), grid; 
hold on; 
end; 
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hold off 

title('Plot of Alpha versus time’) 
xlabel(’Time (sec)'); 
ylabelC Angle of Attack (degrees)'); 
pause;clg; 

% 

axis([0 5 -20 100]); 
for k=l:length(Delta) 

eval([’X_theta=X_theta',num2str(k), ’;']); 
plot(time,X_theta(2,; )),grid; 
hold on; 
end; 
hold off 

title('Plot of Theta versus time') 
xlabel(’Time (sec)’); 
ylabel(’Theta (degrees)'); 
axis(’normal') 
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APPENDIX 2- TERMINAL PHASE MATLAB PROGRAM 



% terminal.m 

% TERMINAL MISSILE-TARGET ENGAGEMENT 
% TIM MULL 
% 24 APRIL 1992 1130 

% This model does not take drag or lift into consideration 
% -uses two different methods to compute sigma derivatives: 

% 1. calculates sig_dot, sig_ddot using analytical expressions 
% 2. uses an observer to estimate sighat,sigdothat,sigddothat 
% 3. controls using sigmadothat,sigmaddothat 

% 4. utilized bang bang on sigmadot,sigmaddot with a small linear region 
% in order to alleviate controller chatter once lock has been achieved. 
clear;clg; 

clc; 

disp([' ’]); 
for k=l:4 

disp(['Control= ’,num2str(k)]); 

control=k; 

dt=0.01; 

(^********** Initial Conditions 5^ jjc 5^ j|c j|c j|c j|c 5|c 5|C 5|C 5|C 
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xm=0; 

vxm=2*340; 

axm=0; 

ym=1000; 

vym=0; 

aym=0; 



% missile horizontal position in meters 
% missile horizontal velocity in m/s 

% missile horizontal acceleration in m/s'^2 
% missile vertical position in meters 
% missile vertical velocity in m/s 

% missile vertical acceleration in m/s‘^2 



xt=5000; 

vxt=-2.5*340; 

axt=0; 

yt=30; 

vyt=0; 

ayt=0; 

time=0; 

- — - 



% target horizontal position in meters 
% target horizontal velocity in m/s 

% target horizontal acceleration in m/s'^2 
% target vertical position in meters 
% target vertical velocity in m/s 

% target vertical acceleration in m/s'^2 



%observer stuff 

pl=-100; 

p2=-100; 

p3=-100; 

K=conv([l -pl],conv([l -p2],[l -p3])); 

kl=K(4); 

k2=K(3); 

k3=K(2); 

X=[0;0;0]; 

% 
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for i=l: 1200 



%loop 1 



%if rem(i,100)==0 
%clg; 

%subplot(221 ),plot(xm,ym) 

%subplot(222),plot(time(l:length(U)),U) 

% end; 

gamma(i)=atan2(vym(i),vxm(i)); % gamma is angle of missile velocity 

% vector with respect to inertial ref 

Vm(i)=sqrt(vxm(i)''2+vym(i)''2); % Vm is missile velocity magnitude 



X(i)=xt(i)-xm(i); 

Y(i)=yt(i)-ym(i); 

VX(i)=vxt(i)-vxm(i); 

VY(i)=vyt(i)-vym(i); 

AX(i)=axt(i)-axm(i) ; 

AY(i)=ayt(i)-aym(i); 

R(i)=sqrt(X(i)'^2+Y(i)'^2); 



% X is relative horizontal position 
% Y is relative vertical position 
% VX is relative horizontal velocity 
% VY is relative vertical velocity 
% AX is relative horizontal acceleration 
% AY is relative vertical acceleration 



Rdot(i)=sqit(VX(i)'^2+VY(i)''2); 

% sigma, sigmadot,sigmaddot are the line 
% of sight between the missile and the 
% target, and its 1st 2 time derivatives 



% This simulation is using analytical expressions for sigma * 
% sigmadot, and sigmaddot. Since the target's velocity and * 
% acceleration cannot be measures directly, a more realistic * 
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% approach is to use an observer to estimate sigmadothat and * 
% sigmaddothat. * 



sigma(i)=atan2( Y (i),X(i)); 

sigmadot(i)=(X(i)*VY(i)-Y(i)*VX(i))/R(i)'^2; 

sigmaddot(i)=(X(i)*AY(i)-Y(i)*AX(i))/R(i)^2... 

-2*(X(i)*VX(i)+Y(i)*VY(i))*(X(i)*VY(i)-Y(i)*VX(i))/R(i)M; 






% 3rd order observer to obtain sigma.dot.hat,sigma.ddot.hat 
A=[-k3 1 0 
-k2 -2*Rdot(i)/R(i) 1 
-kl 0 0]; 

B=[k3 

k2 

kl]; 

[phi ,dd] =c2d( A ,B ,dt) ; 



X(: ,i+l )=phi*X(: ,i)+del*sigma(i); 

sigmadothat(i)=X(2,i); 

sigmaddothat(i)=X(3,i); 






% CONTROL SECTION 



% 



% PROP NAV CONTROL 
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if control==l 



% if loop 1 



U(i)=Vm(i)*4*sigmadot(i); 
if abs(U(i))>15*9.81 % if loop 2 

U(i)=15*9.81*sign(U(i)); 
end; % end if loop 2 

% 

% PROP NAV CONTROL based on sigmadothat 
elseif control==2 % if loop 1 

U(i)=Vm(i)*4*sigmadothat(i); 
if abs(U(i))>15*9.81 % if loop 2 

U(i)=15*9.81*sign(U(i)); 
end; % end if loop 2 

% 

elseif control=3 % if loop 1 

% Bang-Bang control based on sigmadothat, sigmaddothat 

N=15*9.81; % N is maximum normal acceration in Newtons 

test(i)=sigmadothat(i)+(sigmaddothat(i)*abs(sigmaddothat(i)))/(2*N); 

eps=0.001; 

if abs(test(i))>eps % if loop 3 

U(i)=N*sign(test(i)); 

else 

U(i)=N/eps*test(i); 

end; % end if loop 3 

% 

elseif control==4 
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% 



Bang-Bang control based on sigmadot, sigmaddot 
N=1 5*9.81; % N is maximum normal acceration in Newtons 

test(i)=sigmadot(i)-i-(sigmaddot(i)*abs(sigmaddot(i)))/(2*N); 
eps=0.001; 
if abs(test(i))>eps 

U(i)=N*sign(test(i)); 

else 

U(i)=N/eps*test(i); 
end; 
end; 

% 

if rem(i,50)=0 
disp('working'); 
end; 

5fC 5(C 5(C 5(C Sfc SfC 5(C JfC 5fC 5fC HXISSIIO 

% Let Xmsl=[xm(i) State space representation of missile states 

% vxm(i) 

% ym(i) 

% vym(i)]; 

Xmsl(:,i)=[xm(i) 

vxm(i) 

ym(i) 

vym(i)]; 

Am=[0 1 0 0 
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0000 
000 1 
0 0 0 0 ]; 

Bm=[0;-sin(sigma(i));0;cos(sigma(i))]; % control force is perp to sigma!! 

[phim,delm]=c2d(Am,Bm,dt); % convert continuous state-space to discrete 

% update missile states 

Xmsl( : ,i+ 1 )=phim*Xmsl(: ,i)+delm*U(i); 

xm(i-t-l )=Xmsl( 1 ,i+l ); 
vxm(i+ l)=Xmsl(2,i-t-l ); 
axm(i-Hl)=-sin(sigma(i))*U(i); 
ym(i-h 1 )=Xmsl(3 ,i+ 1 ); 
vym(i+l)=Xmsl(4,i+l); 
aym(i-hl)=cos(sigma(i))*U(i); 

3f» j|C 3f» jjc 5|C jJC SjC g ^ ^ jJC SjC jJC 5jC 5jC SjC 5jC 5jC SjC 5j» j|C jjC Sfj 5jC jjC SjC 

% Let Xtgt=[xt(i) 

% vxt(i) 

% yt(i) 

% vyt(i)]; 

Xtgt(:,i)=[xt(i) 

vxt(i) 
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yt(i) 

vyt(i)]; 



At=[0 1 0 0 
0000 
000 1 
0 0 0 0 ] ; 

Bt=[0;axt(i);0;ayt(i)]; 

[phit,delt]=c2d(At,Bt,dt); 

Xtgt(: ,i+ 1 )=phit*Xtgt(: ,i)+delt* 1 ; 

xt(i+l)=Xtgt(l,i+l); 

vxt(i+ 1 )=Xtgt(2,i+ 1 ); 

axt(i+l)=axt(i); 

yt(i+l)=Xtgt(3,i+l); 

vyt(i+l)=Xtgt(4,i+l); 

ayt(i+l)=ayt(i); 

time(i+ 1 )=time(i)+dt; 

%******** check for CPA 

X(i+l)=xt(i+l)-xm(i+l); % X is relative horizontal position 

Y(i+l)=yt(i+l)-ym(i+l); % Y is relative vertical position 
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R(i+l)=sqrt(X(i+l)'^2+Y(i+ir2); 



%******** subroutine to interpolate the 2 closest ranges to find CPA 
ifR(i+l)>R(i) 

xms=[ym(i- 1 :i+l );vym(i- 1 :i+l );xm(i- 1 :i+l ); vxm(i- 1 ; i+1 )] ; 
xts=[y t(i- 1 : i+ 1 ); vy t(i- 1 : i+1 );xt(i- 1 :i+ 1 ); vxt(i- 1 :i+ 1 )] ; 
rl = interp(xms(:,l:2),xts(:,l:2)); 
r2 = interp(xms(:,2:3),xts(;,2:3)); 

Rmin=min([rl r2]); 

break; 

end; 

end; 

j=k; 

eval(['xm',num2str(j),' =xm;']); 
eval([’ym',num2str(j),' =ym;']); 
eval(['xt',num2str(j),’ =xt;']); 
eval(['yt',num2str(j),' =yt;']); 
eval(['Rmin’,num2str(j),’ =Rmin;']); 
eval(['U',num2str(j),’ =U;']); 
eval([’save data’,num2str(j)]); 
end; 

clg; 

clear 

load datal 
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save data5 xml yml xtl ytl 



clear 

load data2 

xm2=xm; 

ym2=ym; 

xt2=xt; 

yt2=yt; 

save dataS xml yml xtl ytl xm2 ym2 xt2 yt2 
clear 

load dataS 

xm3=xm; 

ym3=ym; 

xt3=xt; 

yt3=yt; 

save dataS xml yml xtl ytl xm2 ym2 xt2 yt2 xm3 ym3 xt3 yt3 
clear 

load data4 

xm4=xm; 

ym4=ym; 

xt4=xt; 

yt4=yt; 
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save data5 xml yml xtl ytl xm2 ym2 xt2 yt2 xm3 ym3 xt3 yt3 xm4 ym4 xt4 

yt4 

clear 

% 

% This routine takes data computed in terminal.m and samples it every 
% 8 data point in order to make the graphs less memory intensive 

clear 

load dataS 

xml a=[] ;yml a=[] ;xt 1 a=[] ;y tl a=[] ;xm2a=[] ;ym2a=[] ;xm3a=[] ;ym3 a=[] ; 
xm4a=[] ;ym4a=[] ; 

k=min([length(xml),length(xm2),Iength(xm3),length(xm4)]); 

fori=l:k 



if rem(i,8)=0 

xmla=[xmla xml(i)]; 
ymla=[ymla yml(i)]; 
xtla=[xtla xtl(i)]; 
ytla=[ytla ytl(i)]; 
xm2a=[xm2a xm2(i)]; 
ym2a=[ym2a ym2(i)]; 
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xm3a=[xm3a xm3(i)]; 
ym3a=[ym3a ym3(i)]; 
xm4a=[xm4a xm4(i)]; 
ym4a=[ym4a ym4(i)]; 

end; 

end; 

axis([0 5000 0 1200]); 

plot(xmla,ymla,'-w',xtla,ytla,'-w',xm4a,ym4a,'-w') 

pause; 

plot(xm2a,ym2a,’-w’,xtla,ytla,’-w',xm3a,ym3a,'-w') 

clg;clear 

load datal 
ul=U; 
load data2 
u2=U; 
load data3 
u3=U; 
load data4 
u4=U; 

save control ul u2 u3 u4 
clear 

load control 
dt=0.01; 
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n=inm( [length(u 1 ),length(u2),length(u3),length(u4)] ); 

tfmal=(n-l)*dt; 

time=0.00:dt:tfinal; 

subplot(21 1 ),plot(time,u 1(1 :n),'-w'),title(’control=l ’); 
subplot(212),plot(time,u2(l:n),’-w'),title('control=2'); 
pause;clg; 

subplot(21 1 ),plot(time,u3(l ;n),'-w’),title('control=3'); 
subplot(212),plot(time,u4(l:n),'-w'),title(’control=4'); 
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